#!/bin/bash -e

function info() {
echo Usage: `basename $0` 'in.vcf'
exit 1
}

while getopts  ":p:f:" opt; do
	case  $opt  in
		p) out_prefix=$OPTARG;;
		f) suffix=$OPTARG;;
		*) info;;
	esac
done
shift $(($OPTIND - 1))

if [ $# -lt 1 ]; then info; fi

. $var
# :<<mark
# cat $1 | $java_run/snpsift filter "( GEN[?].AD[?] > 2 )" > $out_prefix.filt_dp.vcf

cat $1 | $java_run/snpsift filter "( GEN[2].DP > 9 )" > $out_prefix.filt_dp.vcf

cat $out_prefix.filt_dp.vcf |$java_run/snpsift filter "( isVariant(GEN[0] ) & ( isVariant(GEN[1]) ) & ( isVariant(GEN[2])) )" > $out_prefix.filt_gt.vcf

java $j_mem -jar $snpsift rmInfo $out_prefix.filt_gt.vcf AF > $out_prefix.rm_af.vcf

java $j_mem -jar $snpsift \
annotate \
-info AF $data_path/gatk/vcf/1000G_phase1.snps.high_confidence.b37.vcf  \
$out_prefix.rm_af.vcf \
> $out_prefix.af.vcf



java $j_mem -jar $snpsift \
annotate \
-info COMMON,G5 $data_path/ncbi/dbsnp/All_20150605.vcf.gz \
$out_prefix.af.vcf \
> $out_prefix.common_g5.vcf

cat $out_prefix.common_g5.vcf | $java_run/snpsift filter "( ! (AF > 0.01) )" > $out_prefix.filt_af.vcf

cat $out_prefix.filt_af.vcf | $java_run/snpsift filter "( ! (COMMON = 1) )" > $out_prefix.filt_common.vcf

cat $out_prefix.filt_common.vcf | $java_run/snpsift filter "( ! (exists G5) )" > $out_prefix.filt_g5.vcf



echo;echo;echo snpeff ann snpeff2gatk
java $tmp -jar $snpeff \
	-c ${snpeff_path}/snpEff.config \
	-v \
	-lof \
	-stats $out_prefix.snpeff_genes.html \
	$snpeff_db_version \
	$out_prefix.filt_g5.vcf \
	> $out_prefix.snpeff.vcf



# java $tmp -jar $snpeff count -i wes.bed $out_prefix.snpeff.vcf > $out_prefix.count.txt



cat $out_prefix.snpeff.vcf |$java_run/snpsift filter "( isHet(GEN[0] ) & ( isHet(GEN[1]) ) & ( isHom(GEN[2])) )" > $out_prefix.filt_gt1.vcf

cat $out_prefix.filt_gt1.vcf | $java_run/snpsift filter "( ANN[0].IMPACT = 'HIGH' )" > $out_prefix.high1.vcf

cat $out_prefix.filt_gt1.vcf | $java_run/snpsift filter "( ANN[0].IMPACT = 'MODERATE' )" > $out_prefix.moderate1.vcf

$java_run/snpsift extractFields $out_prefix.high1.vcf "ANN[0].GENE" |sed '1d'|sort|uniq -c|sed 's/^  *//; s/ /\t/g'|sort -k1nr > $out_prefix.genes.high0.txt

$java_run/snpsift extractFields $out_prefix.moderate1.vcf "ANN[0].GENE" |sed '1d'|sort|uniq -c|sed 's/^  *//; s/ /\t/g'|sort -k1nr > $out_prefix.genes.moderate0.txt


cat $out_prefix.filt_gt1.vcf|$tools/snpEff/scripts/vcfEffOnePerLine.pl > $out_prefix.one_per_line.vcf

cat $out_prefix.one_per_line.vcf | $java_run/snpsift filter "( ANN[0].IMPACT = 'HIGH' )" > $out_prefix.high2.vcf

cat $out_prefix.one_per_line.vcf | $java_run/snpsift filter "( ANN[0].IMPACT = 'MODERATE' )" > $out_prefix.moderate2.vcf

$java_run/snpsift extractFields $out_prefix.high2.vcf "ANN[0].GENE" |sed '1d'|sort|uniq -c|sed 's/^  *//; s/ /\t/g'|sort -k1nr > $out_prefix.genes.high2.txt

$java_run/snpsift extractFields $out_prefix.moderate2.vcf "ANN[0].GENE" |sed '1d'|sort|uniq -c|sed 's/^  *//; s/ /\t/g'|sort -k1nr > $out_prefix.genes.moderate2.txt

# mark

echo;echo;echo snpsift annotate clinvar
java $j_mem -jar $snpsift \
annotate \
$data_path/ncbi/clinvar/clinvar_20150901.vcf.gz \
$out_prefix.moderate2.vcf \
> $out_prefix.clinvar.vcf

echo;echo;echo snpsift annotate gwasCat
java $j_mem -jar $snpsift \
gwasCat \
-db ${data_path}/snpeff/gwascatalog.txt \
$out_prefix.clinvar.vcf \
> $out_prefix.gwas.vcf


. $cmd_done

select_indel_vcf.sh -n -p$out_prefix $out_prefix.filt_g5.vcf
select_indel_vcf.sh -p$out_prefix $out_prefix.filt_g5.vcf

v=$out_prefix.snp.vcf
o=$out_prefix.snp.filt_anno.vcf

java $j_mem -jar /mnt/ilustre/app/medical/tools/GenomeAnalysisTK-3.3-0/GenomeAnalysisTK.jar -R  $ref_genome -T VariantFiltration -o $o --variant $v --filterExpression "FS > 60.0" --filterName "FS_more_than_60" --filterExpression "MQ < 40.0" --filterName "MQ_less_than_40" --filterExpression "HaplotypeScore > 200.0" --filterName "HaplotypeScore_more_than_200"

v=$out_prefix.indel.vcf
o=$out_prefix.indel.filt_anno.vcf

java -Xmx20g -Xms20g -Djava.io.tmpdir=./.tmp -jar /mnt/ilustre/app/medical/tools/GenomeAnalysisTK-3.3-0/GenomeAnalysisTK.jar -R $ref_genome -T VariantFiltration -o $o --variant $v --filterExpression "FS > 200.0" --filterName "FS_more_than_200"

combine_vcfs.sh -p$out_prefix $out_prefix.snp.filt_anno.vcf $out_prefix.indel.filt_anno.vcf



cat $out_prefix.combined.vcf | $java_run/snpsift filter "( FILTER = 'PASS' )" > $out_prefix.pass.vcf

